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ABSTRACT 

Thermal history after inflation is studied in a chaotic inflation model with su- 
persymmetric couplings of the inflaton to matter fields. Time evolution equation is 
solved in a formalism that incorporates both the back reaction of particle production 
and the cosmological expansion. The effect of the parametric resonance gives rise 
to a rapid initial phase of the inflaton decay followed by a slow stage of the Born 
term decay. Thermalization takes place immediately after the first explosive stage 
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for a medium strength of the couphng among created particles. As an apphcation 
we calculate time evolution of the gravitino abundance that is produced by ordinary 
particles directly created from the inflaton decay, which typically results in much 
more enhanced yield than what a naive estimate based on the Born term would 
suggest. 
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I. Introduction 

Inflationary universe models have several attractive features it serves to re- 
solve some outstanding problems in the standard big bang model; how the homo- 
geneous and isotropic universe may emerge from a clumpy spacetime. But how the 
inflationary stage is terminated and smoothly joined to the hot big bang is not fully 
understood at present. In particular, almost all the entropy in our observable part 
of the universe is believed to be created by dissipative oscillation of the inflaton field 
around the minimum of a potential after the rapid cosmic expansion. It is thus 
imperative to deepen our understanding of the dissipative process of infiaton decay 
in order to know more about the matter content of our universe. 

The theory of particle production and associated field dissipation has been much 
advanced recently 0, ||^, [Q, by proper recognition of importance of the 
phenomenon of the parametric resonance. According to the new picture, one must 
distinguish the two processes after infiation; rapid particle production caused by 
oscillation of the coherent infiaton field and subsequent thermalization process. The 
thermalization occurs due to interaction and decay of produced particles, which 
involves some weak coupling constant. In contrast, an enormously large number 
of particles is created prior to particle interaction by the oscillating field due to 
the instability in the parametric resonance region, if a combination of the initial 
oscillation amplitude of the infiaton and its coupling to bose matter fields exceeds a 
critical value. The time scale of resonant particle production and associated infiaton 
decay is typically of the order of the inverse infiaton mass, which is short. Thus the 
old picture based on the Born term decay of infiaton must be modified. 

In this paper we investigate thermal history after infiation and demonstrate that 
Planck scale remnants such as the gravitino, having rates ordinarily suppressed by 
the Planck mass, may be overproduced from those directly produced by the infiaton 
decay. It has been argued ||^ for some time that if infiation ends with a low enough 
reheating temperature, the gravitino overproduction may be avoided in order not 
to upset the success of the standard big bang cosmology. We compute here the 
gravitino abundance right after infiation, starting from the first principle of particle 
physics. It is found that in a wide range of parameters copious production of ordinary 
particles in an early epoch of the infiaton decay gives rise to a significant amount 
of gravitino abundance and subsequent dilution at later epochs is not very effective. 
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Our result on the gravitino abundance may equally be applied to any exotic particles 
that have coupling to ordinary matter with the gravitational strength. Thus one may 
entertain a possibility that stable secondary remnants created by direct products of 
the inflaton decay, being capable of probing the Planck scale physics, may actually 
comprise the bulk of the present mass density of the universe as the non-baryonic 
dark matter. If the gravitino is the Lightest Supersymmetric Particle (LSP), there 
is a window in the parameter space for which the gravitino dominated universe gives 
the solution to the nagging gravitino problem. 

Although the complete theory of reheating still does not exist, it is necessary to 
extend previous works on the inflaton decay and associated particle production in 
order to discuss relic abundance of the Planck scale remnant. A fundamental picture 
of the new theory of reheating was outlined by Kofman, Linde, and Starobinsky 
0]. In our previous papers 0, we discussed in great detail a field theoretical 
framework of the inflaton decay and worked out time evolution, taking into account 
the back reaction of produced particles. In the following we first explain need to 
study in more detail both the quartic and the Yukawa coupling cases in a class 
of supersymmetric models. We then give our own account of the inflaton decay, 
highlighted by a semi-classical set of differential equations for the inflaton energy 
and the created radiation energy density in the expanding universe. 

Our present work, when limited to the inflaton decay, improves over previous 
works mathematical treatment of the parametric resonance in the large amplitude 
limit and incorporates the small amplitude region at the same time. Moreover, by 
integrating over short time behavior we derive tractable time evolution equation 
for the inflaton and the created radiation energy density, from which some physical 
consequences of the inflaton decay, otherwise ^ difficult to derive, can be extracted. 
Both the back reaction and the cosmological expansion play crucial roles in our 
analysis. Physical consequences that come out of our analysis are not much different 
at the qualitative level from those of Ref in the first stage of infiaton decay, 
but differ substantially in later stages. Most notably, we find unlike others that 
thermalization takes place immediately after the first explosive stage of the infiaton 
decay, both because in our models created particles interact among themselves with 
ordinary strength of the medium magnitude, and also because tremendous increase 
of the particle number density occurs abruptly at a time in the first stage. 

Relic gravitino abundance from infiation measured by the ratio of the number 
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densities, n2,/2/n^ is roughly given by y^/'mpi with the number density of 
produced massless particles and T a typical particle energy. This ratio agrees with 
the usual formula |P] Tfj/mpi, valid when thermalization is assumed with a reheat 
temperature Tr. In the present work we first check how and when thermalization 
condition is met and then follow time evolution of the gravitino abundance after 
thermalization takes place. Fortunately the maximal particle production and ther- 
malization occurs almost at the same time in an early epoch of the inflaton decay 
for a wide range of parameters, which makes our analysis easier. By solving time 
evolution equation, one does not have to separate epochs of gravitino production 
and its subsequent dilution, sometimes difficult to distinguish. It turns out that 
the gravitino abundance thus computed cannot be given by a simple formula like 
Tfl/mpi. 

It would be easier to understand the rest of this paper if we first sketch essence of 
time evolution that comes out of our subsequent detailed analysis. There are three 
important stages of the infiaton decay. After an initial period of the infiaton domi- 
nance, abrupt decrease of the infiaton amplitude occurs, accompanying an enormous 
and the maximal rate of particle production. Due to a high rate of interaction among 
produced particles thermalization takes place immediately after this time. The sec- 
ond stage is a slow process of particle production, that continues until the last stage 
of complete decay when the Hubble time becomes comparable to the infiaton decay 
time given by the Born formula. Difference in the first two stages lies in the relevant 
6 value (precisely defined later by a dimensionless coupling in the mode equation of 
quantum field to be created), large 6 in the first stage of rapid decay and small in 
the second stage of slow decay. Although the bulk of the universe may be taken ther- 
malized after the first explosive stage, it would be misleading to assign a single reheat 
temperature. Later we shall define two temperatures, the initial reheat temperature 
Tj immediately after the the explosive stage, and the final reheat temperature Ty 
after the complete Born decay. For a large enough coupling the Planck scale relic is 
most abundantly produced in the aftermath of the first explosive decay, followed in 
the second stage by a slow dilution of the abundance until the Born decay. But this 
picture is modifed both for small coupling and for extremely large coupling. 
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II. Model 



Let C,{t) a generic inflaton field treated here as a classical and homogeneous vari- 
able. Coupling of the inflaton to a generic bose matter field if is assumed to contain 
both the Yukawa (^v^^) and the quartic (^^v^^) type interaction, as so happens in 
supersymmetric models we consider below. Under a periodic oscillation the quan- 
tum state of if system exhibits instability due to the presence of the exponentially 
growing mode in the parametric resonance region. Evolution of the quantum ground 
state is simple 0, 0: the quantum state is described by a density matrix that 
has rapidly oscillating off-diagonal elements. Hence the averaged density matrix 
smoothed over short time scale takes a nearly diagonal form. Thus the system is 
almost classical: even a small bit of interaction among created particles is expected 
to break the coherence and yield the classical behavior of a mixed state. Instability 
of the quantum if system is thus interpreted to give rise to particle production and 
associated inflaton decay. This occurs in a large region of the parameter space of 
the coupling constant and the initial ^ oscillation amplitude. We use these results of 
particle production from the ground state even to the case from mixed states such as 
the thermal state at late stages. The basic attitude we take here is that a real pro- 
cess of particle production may be insensitive to the state in which it occurs, while 
quantum fluctuation involving virtual states of all kinds may sensitively depend on 
the state under consideration. 

For definiteness we concentrate in this work on the simplest kind of inflation 
model, the chaotic inflation [|T^ with a mass term potential |m|^^. It would be 
instructive to give beforehand a crude range of parameters we have in mind. The 
fluctuation amplitude in the density contrast observed by COBE suggests the inflaton 
mass in the range of 10^'^ GeV [11|, while the initial oscillation amplitude 
is of order the Planck scale. Whenever we need concrete numbers, we use = 
-y/^mpi ^ 6 X lO^^GeV and ^q/tti^ = 5 x 10^. The quartic coupling constant g"^ and 
a dimensionless combination of the amplitude 9q ^ g'^^o/fn^ are the key parameters 
for subsequent evolution after inflation. For g"^ ^ ^|/^o ~ 10~^^ the large amplitude 
oscillation with ^ 1 is realized, thus our analysis in the following is applicable in 
a large region of coupling constant. 

A class of supersymmetric models we consider is characterized by the superpoten- 
tial of a generic form, W = ^ ^^^"^ ~'~ ^ 5'^'^^ ' ^i^^ 'P representing a chiral multiplet 
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and ^ the inflaton field. The magnitude of the coupling constant g is not much 
constrained and can be very small. This gives ordinary potential terms of the form, 
I 9^ i'^v'^ + \ 9^i^f'^ ; foi' the ^ coupling. Note that both the Yukawa and the quartic 
coupling terms are here with a definite relation between the two. The number of real 
scalar fields that couple to the inflaton field ^ is taken to be 8 (the number of real 
components of the two Higgs doublets) as in the minimal supersymmetric model. 
Although we do not write interaction terms of bose fields ip among themselves, they 
are also induced by other terms of the Lagrangian besides the one derived from this 
superpotential. We assume that these generic bose fields have interaction strength 
of medium size, 5'^/(47r) = 1 — 10~^. These interactions are important to establish 
rapid thermalization once copious particle production occurs. A quantity to mea- 
sure strength of the resonant particle production is the variable mass term in the 
mode equation for the (p field. Ratio of the variable mass terms derived from the 
quartic and the Yukawa coupling is ~ Q^/^i-, which is large in the initial stage of 
the inflaton decay. Thus with an initial amplitude of g^o/m^ g ■ 10^ ^ 1 , it is 
sufficient to consider the quartic coupling alone initially, although the Yukawa cou- 
pling is important in later stages and ultimately responsible for the Born term decay 
of ^ — > <y9 . It would be worthwhile to recall the often cited reheat temperature 
of O[0.1 ] mpiF^ ^ gVNd ■ 10^^ GeV with a number of decay channels, in the 
chaotic inflation model with = 10^^ GeV . The problem with this naive formula 
is that it is not clear at what stage of the inflaton decay one should apply this reheat 
temperature, especially because the new theory of inflaton decay implies a prolonged 
period of particle production lasting until the complete Born decay. 

III. Parametric Resonance and Particle Production 

Analysis of the instability is divided into the small 6 and the large 6 regions, 
basic formulation of which was respectively given in Ref |^ and in Ref . With the 
effective coupling term 2g'^^'^ip'^ specified in the initial stage, the Fourier mode 
of the wave vector k obeys 'fk +{k'^ + 4(yf^^^)v^fc = 0, aside from effects of the 
cosmological expansion. The effective ^^v'^ vertex is different from the bare vertex 
derivable from the bare Lagrangian, due to additional two diagrams that contribute 
in the low energy limit with equal strength to the bare vertex. Taking the sinusoidal 
oscillation of the ^ inflaton field into account leads to the standard form of Mathieu 
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type equation, 

<Pu 

— + (/i-2^cos(2r)) M = 0, (1) 
n = — ^ H ^ , 6 = — 5- , r = md , (2) 

with ^ the amphtude of oscillation. The relevant parameter region is thus restricted 
to h — 26 > for the quartic coupling. In contrast the Yukawa coupling term 
becomes dominant in the late stages and gives 

h = —,e = , T = -m^t, (3) 

m| 2 

thus the relevant region is not restricted to h — 26 > 0. The parametric resonance 
occurs in an infinite number of bands of the instability region separated by stability 
bands. In the instability band the mode function exhibits the exponential growth, 
e^'^ with time r, and accordingly ip particle production proceeds. Thus important 
quantities one must extract from mathematical analysis of the mode equation are 
the growth rate A and the phase space volume that contributes to the resonant ^ 
decay. The rate A in general depends on the energy u of created particle, but for a 
global analysis we perform in this work, it is sufficient in some cases to know a value 
of A averaged over the relevant phase space. 

As was shown elsewhere 0, in the large amplitude limit of ^ 1 one can 
replace the sinusoidal function in the potential term of the Mathieu equation by a 
periodic potential, piecewise quadratic. The exponential growth rate X of ip particle 
production in the resonance region is then given by 

A = - ln(v^+ V;^^^), X = (1 + e-^^^) cos^V', (4) 

r ((1 -iv^e)/2 



^lJ = —V9 + V9e\n{n9^/*) + - $5 In 



r((l + ^v^e)/2 



(5) 



with e = ^ — 1 . Even for a 6' as small as 0[1] the formula above is a qualitatively good 
approximation to the growth rate. The parametric resonance region is determined 
by a; > 1, which is separated by the stability gap of x < 1. For instance, one may 
plot the growth rate A given by Eq (^) against e for a fixed 6, as illustrated in Fig.l. 
Behavior of the growth rate A is very much different, depending on whether h > 29 01 
h < 26, in other words, e > or e < 0. The relevant coupling in the large 6 region in 
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our supersymmetric model is the quartic coupling, which necessarily has e > 0, hence 
we shall concentrate on this region here. As a positive e increases, the peak growth 
rate A^ax within an island of the instability region decreases, but the first peak A 
value always remains of order 0.2. This is related to that at e = 0, or along h = 26, 
the maximum value analytically calculated is given by Amax = ln(v^+ l)/7r ~ 0.28 . 
We numerically computed the average value of the growth rate in the e > region, 
taking a wide range of 6, 10 — 10^. The average of the first A peak value is ~ 0.15. 
On the other hand, the phase space volume element multiplied by energy to get the 
ip energy density is 2g^^'^/TT'^ times (0.07 — 0.1)/6' numerically. 

With cosmological evolution the 6 value varies, and at late stages of the inflaton 
decay one must consider the small 6 region, even if it starts in the large 6 region 
initially. Actually there are two causes of 6 variation: either by decrease of the 
oscillation amplitude ^ or by the C, mass {m^) variation, which we later turn to. 

The most important resonance band in the small 6 region, and the only one we 
consider in this paper is the first band which is approximately given, with 6 < 1, by 

i 

Max[0,h^]<h<h+, h± = l±e , (6) 

8 

with h and 6 related to physical quantities of the Yukawa coupling. The phase 
space of the first band as measured by the h ( oc cj^) width becomes nonvanishing at 
6 ^ \/24 — 4 ~ 0.90, diminished at Vl52 — 12 ~ 0.33 , until it becomes very small 
like 26 for very small 6. 

There is an important limitation in applying the resonant decay when cosmolog- 
ical expansion is taken into account: the cosmological redshift may push relevant 
frequency out of the resonance band 0, 0]. The effect is particularly important in 
the second stage of slow decay. The growth factor is roughly modified as follows. 
Time interval during which a particular mode stays within the first band at the 
Hubble time 1/H is given by At = e/{2H) = 6/H , using the width of the first 
band e. Multiplied by time averaged rate A ~ 6/A, this gives the growth factor in 
the exponent, ^ ■m^6'^ / {8H) . More precise estimate is possible by integrating the 
growth rate while remaining within the narrow band during the frequency change 
Au; ~ uHAt, 

^ fdtX = ^ fdAu ^^h-~ fh-]tl = ^^^ (7) 
2 J A J 16 if ' ^ ' 

I at I 
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where Acj = uj — tUj-es with cijj.es = the center position of the resonance. The 

relevant phase space is estimated similarly. For instance, it contributes to the created 
energy density, 

1 o ,(iAcj, mf , ^ 

■ujU^ = —^H. 8 



27r2 '"'^ dt ' 327r2 
Thus with Ref we may use for the source term of radiation. 



m 



4 ^ ^ q2 



TT 171^9 



,i7sinh^[ ^1. (9) 

327r2 ^ 16 ^ ^ ^ 

There is a problem with using this formula all the way down to the smallest value 
of the amplitude ^. In the small 6 = 2y/2 g^yn^/m^^'^ limit the number density of 
^ particles is very small, and the coherence effect is not important. One then expects 
that the Born rate for a single bosonic decay channel, = g^m^^'^ /{GAn) , should 
appear as the source term, leading to the ordinary exponential decay law. Indeed, 
it was shown |Q in the case without cosmological expansion taken into account that 
the exponential decay law follows. The source term above however gives a factor 
oc g^ in the small g limit, which cannot be reconciled with the Born rate oc g^. 
We suspect that a more refined treatment of the mode sum for the created particle 
number in the expanding universe gives both the parametric resonance effect in the 
medium small amplitude region and the Born rate behavior in the extreme limit of 
small amplitude. We have not succeeded in deriving this from a single formula in a 
convincing way. In the present work we are content to take a more phenomenological 
point and simply add the two terms, the formula above Eq (^) and the Born rate, 
anticipating that the two relevant regions do not overlap in the parameter space 
and when one effect is large, the other is negligible. This should be a practically 
good approximation to physical processes involved, although a rigorous derivation is 
welcome. 

Dissipation due to particle production necessarily accompanies quantum fluctua- 
tion of the (fi field in the state. As noted in Ref 0, [Q, this has an important effect on 
variation of the ^ mass via g'^{f'^) ■ The expectation value here may be evaluated by 
using an exact relation in each mode between the produced particle number and the 
fluctuation, which follows since the wave function of the quantum state is given by a 
Gaussian form. With the exponentially small term neglected, this relation reads as 
{]^uj) ~ — — - {ill) , with the coordinate variable of harmonic oscillator of mode 
u. In both the large 6 or the small 6 region one can neglect A dependent term here. 
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First, in the large amplitude region typical resonance energy uP' — 0\g^m^ ^ m| . 
Second, in the small amplitude limit Am^ ~ dm^jAL , which becomes smaller than 
itself. In the following we denote the variable mass by m^\x. The mass variation in 
the region where the Yukawa coupling dominates involves one loop two point func- 
tion in which two external ^ lines arc separated unlike the previous tadpole term. 
Computation of this particular diagram in the quantum (/? state in question has not 
been worked out. We replace the entire one loop contribution by the single tadpole 
term (99^) we calculate here, hoping that the neglected one loop term gives a con- 
tribution comparable to, and not too much larger than, the calculated tadpole. If 
this is the case, our neglect of the other contribution would be excused, because the 
mass variation has been found not to be drastic, /i < O[100] in the coupling region 
of our interest. Even more precise mass variation term would not change much of 
our conclusion. 

Summarizing the particle production rate of species (= 8 in the minimal 
supersymmetric model) and the rate of mass variation, we have obtained the mode- 
summed expression: with dimensionless functions defined by 

y = Pi/pl ' z = pr/pl, pI = \ ^Uo > P = m^{t)/m^ , (10) 
j^T.{P.)^l^ieoBsy, g'j^^{^')^2mlAs,, (11) 

in the large 9 region 

, 0.17x0.3 o 6 „ 0.1 X 0.3 2 / X 

= exp[0.3 / dr jji{T')], Su = v^exp[0.3 / dT'jjL{T')], (13) 
JO Jo 

by taking an average value 0.3 for 2A, while in the small 9 region 



Ar^^ Br^ i^f (14) 



y 1 I — \ — • u2r /I y 



Sy^ \x \/y + z sinh [7r6'o . . — ■ — ] , = P\/y + z sinh [tt9o , — ■ — ] , (15) 

H y/y + z li yjy + z 

with understanding that the Born term should be added in the source term of evo- 
lution equation. 
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IV. Inflaton Decay and Thermal History after Inflation 

Time evolution of the inflaton and the radiation energy density is given by a 
coupled set of differential equations. This equation is derived by considering energy 
balance between ^ and (p with the cosmological expansion taken into account. Vari- 
ation of the ^ mass involves the expectation value (y?^) which varies like ^^^(t) , 
as the cosmological scale factor a(r) increases. Coupled evolution of the ^ energy, 
the (/? energy density, and the ^ mass is then described in terms of the dimensionless 
quantities by 

-Bsy-d'yy, (16) 

Bsy + d'yy, (17) 

-)=As,, (18) 
/i 

with d = ^ , and 7 = IN^g^ j (327r) the dimensionless Born decay rate. 
The number of ^ decay channels is doubled here since fermions of supersymmetric 
partners contribute to the Born rate, but not to the parametric resonance effect. 
The initial time for integration may be set at the Hubble H = m^, since the damped 
oscillation effectively starts there. This amounts to letting d = 1, assuming no initial 
radiation present; y{0) = 1 , z{0) = . 

The matching from the large to the small amplitude region is somewhat delicate. 
Recall the deflnition of 6 for the Yukawa and the quartic couplings, which leads to 
a relation, 6y = 2\^/fi. Suppose that the mass variation is small, ~ 1. It is 
then possible to have a situation that despite of 64 < 0.33 the Yukawa 6y > 0.33. 
We should ensure for the matching that the larger one of 6y or 64 is small, and 
within the small amplitude region. Without embarking on a too much elaborate 
analysis, perhaps not worthy of its effort at the present level of investigation, we 
practically adopt as the small amplitude condition Oy < 0.33 and use the large 
amplitude formula for the quartic coupling. For completeness we did also integrate 
the differential equation by a more elaborate matching; in the intermediate region 
with 6*4 < 0.33 and Oy > 0.33 we have replaced the source term by formulas applicable 
in the large Yukawa coupling region, as worked out in Ref . Difference in this case 
is appreciable only in a narrow region of the coupling g, as seen in Fig 7b for the flnal 
gravitino abundance. Hence the flnal outcome does not seem to sensitively depend 
on the matching condition. 



y' + 3dy/y + zy = 
z' + ^d^Jy^^ z = 
n' + d\Jy + z (/i - 
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We did some extensive numerical investigation of time evolution equation for the 
set of parameters, g and ^o/m^. Here we restrict to the special case of ^o/^$ = 5x10^ 
and g — 0[1 — 10~^], consistent with COBE observation in the chaotic inflationary 
scenario. With a coupling less than 10~^ the effect of the parametric resonance is not 
significant. This critical strength for the parametric resonance may be determined 
by the condition, 64 ^ 1, leading to g ^ ^g/'Co ~ 10^^. Moreover, it has been 
found by the present work that not only for g < 10^^ but also for g > 1, the naive 
estimate of the gravitino abundance by the Born formula happens to be numerically 
close to the correct value. This is the reason we concentrate on the coupling region 
of 1 > g > 10^^ . Before we present result of numerical computation, we discuss 
question of thermalization, coherence and how to compute the gravitino abundance. 

Although the process towards complete thermalizaion is not fully understood at 
present, the following argument suggests that with pr ~ at the exit t = tg from 
the first explosive stage, thermalization takes place immediately. For a wide range 
of the coupling, the condition pr ^ p^ holds, for instance for > 3 x 10~^, 10 > 
Pr/p^ > 0.01. Note first that a typical energy of created particles at the exit time, 
not necessarily thermalized as yet, is of order T ^ gm^ p^ ~ 9^^'^ p ~ m^p . 
Interaction rate among (p particles compared with the Hubble rate is very crudely 
estimated as 



with r = Pr/Pc ratio of the radiation to the ^ energy density at the exit. Thus 
with ^0 close to the Planck scale, T ^ H at te since pe < O[10^] in the parameter 
range of our interest, and thermalization takes place almost at the same time as the 
explosive resonant particle production. We did a consistency check of the condition 
of thermalization F/H ^ 1 in our numerical analysis, by using the interaction rate 
r = aT with a > 10~^ and T the temperature estimated from the radiation energy 
density Pr- It turns out that this thermalization condition is always obeyed after the 
explosive exit from the large amplitude region, as shown in Fig 2. 

Even the multi-particle reaction that works very efficiently as the energy exchange 
process may occur frequently. For instance, rate of the (N+l) — > 2 body reaction 
behaves on dimensional grounds like 




(19) 



) 



N 



(20) 
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with as — g1/{i.T:) a typical medium size coupling constant, n the number density 
of Ng species of colliding particles, T an average energy of colliding particles. The 
point of this discussion is that in situations immediately after the explosive stage, 
but before thermalization, there are more low-energy particles than what is expected 
from the thermal distribution, which implies that n ^ T^. With ^ 3> 1, one 
may expect that both the chemical and the thermal equilibrium is established rapidly. 
The number of reactions needed to reach equilibrium may crudely be estimated as 
follows. A typical amount of exchanged energy per reaction is very roughly T m^/Xe 
as noted above. To get the thermal energy of order p^^, one needs energy exchange 
reactions Nj. times, with 



under the assumption of ~ at the exit. This is not a huge number, at most 
100, and one may expect that equilibrium is reached fast. 

Time evolution of the thermal temperature T, defined from the radiation energy 
density by T = [pr/{^N)Y/'^ with N — 200, is shown for a few choices of the 
coupling in Fig 3, together with evolution of the net entropy sa^{t)] the true mea- 
sure of dissipation. In all these cases plateau regions of the net entropy are clearly 
observed indicating period of quiet evolution, but behavior of the transitional era is 
different depending on the coupling g. Although the temperature is calculated from 
the radiation energy density at all times, thermalization does not take place prior to 
the explosive stage, hence the temperature in this early epoch has no meaning. 

Despite of our detailed control of time evolution it would be useful to present 
what is left in the end after the second reheating due to the Born decay. The final 
reheat temperature defined at a time when the net entropy pl^^ settles down to 
its asymptotic value (effectively replaced by that at a time, 10 x Born lifetime), is 
thus plotted against the coupling g in Fig 4. Remarkably, this is not much different 



from the usually estimated reheat temperature of order 0.1 y mpi , using the Born 
formula. We however warn that this does not necessarily support the simple picture 
of the inflaton decay based on the Born term: in the intermediate times between 
the explosive stage and the Born stage there exits a prolonged epoch with higher 
temperatures. Physical processes that involve higher temperatures such as baryo- 
genesis or gravitino production must be studied taking into account this period. For 
reference we plotted in Fig 5 the initial reheat temperature Tj defined immediately 




(21) 
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after the explosive stage. The maximal temperature we observe in the whole range 
of the coupling g is ^ 90m^ ^ 10^^ GeV close to the GUT scale, right after the 
explosive stage. As repeatedly stressed here and will be explicitly demonstrated 
later for the gravitino abundance, the correct thermal history cannot be described 
in terms of a single temperature, either Tj or Tj. 

Despite of immediate thermalization after the explosive first stage decay, energy 
balance equation derived without considering effect of thermalization is presumably 
little affected since the energy conservation holds globally irrespective of energy 
exchange process towards thermalization. On the other hand, the expectation value 
{(f'^) may differ considerably in thermal equilibrium from that we derived previously, 
which may then modify the equation for the mass variation. For simplicity we ignore 
this thermal effect in the present work. 

Coherence of the cosmic ^ field is crucial to the phenomenon of the parametric 
resonance, especially in the large amplitude region. One may estimate in the fol- 
lowing way how and in what parameter region the coherence is lost. First, one may 
regard the homogeneous and coherent field as an aggregate of massive particles at 
rest with a phase coherence precisely tuned. Suppose that some of these ^ particles 
interact with created (p particles, which partially breaks the precise phase coherence. 
Thus whenever interaction rate between ^ and (p particles dominates over the 
cosmic expansion rate H, the coherence is considered to be lost. To estimate the 
interaction rate F^ = (av) ritpN with N — 200 , one first computes the cross section 
a. As to the relevant energy of (p particles the thermal average, ~ 2.77", may be 
taken, which becomes as large as ~ 100 in some situations we consider. Thus for 
simplicity we take the high energy limit cu ^ of the cross section. This results in 

4 

the cross section, — . The ^ interaction rate relative to the Hubble rate is then 
estimated as 

^ ^ 0.68^ VT5N g"^^ 1.2 X 10^ g\ (22) 
H 27r* 

assuming the radiation dominated thermal bath of temperature T. Thus one has to 
worry about the loss of coherence only for a large coupling of order 0.05 or larger. 
We have checked more carefully at what coupling g the coherence is lost, by using 
numerically computed Hubble rate. In most of the parameter region except for 
g > 0.08 the ratio of the ^ interaction rate to the Hubble rate does not exceed unity. 
In our presentation of various output quantities we ignore the coherence constraint. 
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since for these large couplings of g > 0.08 the gravitino abundance itself is too large 
even if the parametric resonance effect is altogether neglected. 



V. Gravitino Abundance 



Time evolution of the gravitino abundance follows the rate equation, 

dn 



-^ + 3Hny2 = {^v)nl, (23) 

since the annihilation term may be neglected due to the low rate of production, as 
consistently checked by the outcome, <^ . The cross section S contains all 
contributing channels, and it was computed as S = 250/mpi |]12[. The remaining 
uncertainty of the (p number density n<^ is fixed by using the relation to the energy 
density in a thermal bath, valid at times after the explosive stage; 

03) 30p, 3/4 

We plotted in Fig 6 time evolution of the gravitino abundance for a few choices of the 
coupling. An unexpected result is that for large enough couplings a sizable gravitino 
yield produced immediately after the explosive stage is almost completely canceled 
by the late Born decay, resulting in the final abundance nearly equal to the naive 
rate estimated by the Born term. 

The final gravitino abundance against the coupling is shown in Fig 7, together 
with the abundance obtained without considering the resonance effect. In other 
words, we also evaluate for comparison time evolution of the gravitino abundance 
with the Born decay term 7^ as the source, without assuming either the instan- 
taneous decay or the instantaneous reheating. This way of computing the gravitino 
abundance would give more precise prediction of the Born term decay than the usual 
simple estimate based on the instantaneous decay and reheating. The final gravitino 
abundance is substantially larger in the range of couplings, 10~^ < 9 < 0.1 than what 
the Born decay would suggest (approximately 8 x 10^^ g). It should be noted that 
gravitino production and its subsequent dilution is a process that goes on continu- 
ously from the first explosive stage to the final stage of the inflaton decay. Although 
the final reheat temperature is not too much different from the Born estimate, what 
happens in the intermediate stages may differ for different couplings, and the final 
abundance has a complicated behavior as a function of the coupling. 
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Fate of the produced gravitinos has been much discussed in the hterature. There 
are strong constraints from the decay products of the photon and LSP (Lightest 
Supersymmetirc Particle) |T^. Production of both the electromagnetic and the 
quark showers [|r^ around the epoch of nucleosynthesis must be much suppressed in 
order not to destroy produced light elements. The decay mode containing showers 
of strongly interacting particles is very much model dependent due to the unknown 
scalar quark and the gluino mass spectrum, while photo-dissociation of D, ^He, ^He 
gives a strong constraint on the abundance s with s the entropy density if the 



gravitino mass is less than a few TeV ||12 



Here we merely mention the constraint from a possible overclosure of the 
present mass density of the universe due to the stable LSP denoted here by Xi 
coming from the gravitino decay. Recall first the decay rate, valid for the gravitino 
mass m3/2 <^ M = m^i/ -\/87r, 

T{^, ^ + A,) = ^ ^ ~ (4 X 10« sec)-^ (m3/2/100 GeV)^ , (25) 

where the specified decay mode is the gauge boson {Ai,) and the gaugino (x)- Relic 
LSP abundance n^/ s is given either by the abundance of the parent gravitino s 
if this yield is not too large, or by the thermal annihilation abundance if the initial 
abundance is large, 

!^ = Min (^,(!^)_) , (^)_~3.8,;V2^ (26) 

S \ S S J S {(y^nnV) Td 

The annihilation abundance (n^/s)ann has been estimated from the balance between 
the annihilation and the Hubble rates at the freeze-out T = To, ((Xann v) = H(Td) , 
with the number of massless species at the freeze-out. 

Let us take as an example the LSP mass = 100 GeV. From Fig 7 we find that 
the condition of non-overclosure with = 11^/2 implies g < {1 — 2) x 10~^. If this 
bound is violated, one must reduce the LSP abundance via thermal annihilation that 
may later occur, in order not to contradict with the closure density. For simplicity 
let us take a likely possibility of x being the bino, in which case the annihilation rate 
at temperature T may be given by 

T T 

(o-annt') = a + b— + 0{ — f , & = E ^// ~ ^O"' , 2 / 2,2 ' (28) 
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assuming a common mass rrijoi exchanged particles. The condition of non-overclosure 
then leads to 

7713/2 > 100 TeV . (29) 



We took over the result of a similar analysis of Ref [|T6[ , applied for the Polonyi field 
decay. Thus with = O[100GeV], one must have either a very small coupling 
g < 10~^ or a gravitino mass larger than 100 TeV, which is not favored from the 
point of low energy SUSY breaking. The constraint becomes even more stringent 
if one also considers destruction of light elements by electromagnetic and hadronic 
showers due to the gravitino decay. 



More interesting is the possibility of the gravitino being LSP |[T^. The old over- 
closure bound < IkeV |]T7| is based on the thermal abundance of gravitino 
production, hence is not applied at the face value if inflation dilutes away thermal 
gravitinos, as in our case. Besides the primordial production immediately after in- 
flation computed by us, there may be a new source of the gravitino production due 
to the next-to-the-lightest superparticle (NSP) decay. The production rate of this 
process has been found negligible and the corresponding bound on the gravitino 
mass, including argument of photo-dissociation due to radiative decay component. 



have been given in Ref [0. We should however derive the bound, not in terms of 
the reheat temperature Tr and the gravitino mass, but in terms of the coupling g 
and the mass, since the reheat temperature is not well defined in the new theory of 
reheating. The result is shown in Fig 8. In this argument it is important [jl8| to 
modify the gravitino production cross section to 

2 

4.5 , % , (30) 

^3/2 '^pl 

assuming that gravitino production proceeds mainly from gluino, as is the case unless 
the gluino mass rrig is too small. The reason the gravitino mass ^3/2 appears in the 
denominator is that with a light gravitino the production mainly proceeds via longi- 
tudinal component, Goldstino. Although not shown in Fig 8, the closure condition 
^3/2^^ = 1 is met roughly with (7/(7713/2/1 GeV) = 3xl0~^ evenform3/2 < 10~^ GeV, 
down to 7713/2 ~ IkeV at which the thermal, hence the maximal yield is reached. 
Thus with the gravitino mass 777.3/2 < 1 keV, thermal production becomes possible 
even after inflation. Also is not shown the constraint from the photo-dissociation of 
light elements which becomes very stringent for a gravitino mass less than a few GeV 
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121. The maximal reheat temperature right after the explosive stage can become 



as large as ^ 2 x 10^^ GeV, with g ^ Q x 10"^ , 7713/2 ~ 1 GeV (for nig = 150 GeV) 
in the gravitino dominated universe. This high temperature might help to generate 
the baryon asymmetry at the GUT scale. The gravitino dominated universe thus re- 
mains as an exciting solution to the nagging gravitino problem and deserves further 
study. We note that our calculation is presumably the first instance that computed 
the gravitino abundance starting from the first principle. 



In summary, the thermal history after inflation has been elucidated in a simple 
model of chaotic inflation that contains a generic type of coupling to matter fields. 
The outcome of the inflation cannot be described by a single reheating temperature, 
since the initial phase dominated by the parametric resonance effect, ends with much 
higher temperatures followed by a slow inflaton decay due to the usual Born term. 
High energy processes are much affected in intermediate stages, and we explicitly 
demonstrated this by following time evolution of the gravitino abundance. 
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Figure caption 



Fig 1 

Growth rate in the instability region plotted against e = ^ — 1 for a fixed value 
of e = 1000. 

Fig 2 

Thermalization condition measured by the quantity, the thermal temperature 
relative to the Hubble rate T/H that closely approximates the interaction rate, 

Fig 3 

Time evolution of the temperature in the unit of the inflaton mass (solid line) 
and the net entropy defined by the entropy density times the scale factor to the 3rd 
power (dash-dotted) for three cases of the coupling (A) g = 10~^, (B) g — 10~^ (C) 
g — The broken line shows the temperature as computed from the source term 

of the Born decay without effect of the parametric resonance. The time scale is given 
in the unit of the Born decay lifetime t^^—- 

Fig 4 

Final temperature in the unit of plotted against the coupling g. The final time 
is taken at 10 xtb (Born hfetime). Result based on the Born term decay depicted 
by the broken line, coincides with the more precise result (solid line) within precision 
of this figure. 

Fig 5 

Initial temperature in the unit of immediately after the explosive stage, at 
the time of Oy — 0.33 . 

Fig 6 

Time evolution of the gravitino abundance (solid line) plotted against time the 
unit of the Born decay lifetime ^e^^, for three couplings (A) g = 10~^, (B) g = 10~^ 
(C) g = 10~^. For comparison result due to the Born term alone is shown by the 
broken lines. 
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Fig 7A 

Final gravitino abundance relative to the entropy density plotted against the 
coupling g (solid line). The final time is taken at 10 xtb (Born lifetime). For 
comparison is shown the same abundance calculated without the source term due to 

resonance effects (broken line). 

Fig 7B 

The same gravitino abundance as in Fig 7B, with the two step matching as 
described in the text (sohd fine) , to be compared with the one step matching in Fig 
7A (dash-dotted fine) , and with the Born rate (broken fine) . 

Fig 8 

Region of the parameters g,m3/2 in which the stable gravitino overdoses the 
present mass density of universe (shaded region) . The solid lines indicate the closure 
contours with the two choices of the gluino mass rrig, and the broken lines the similar 
ones obtained from the Born source term. 
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